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Abstract 

Sensory neurons are often described in terms of a receptive field, that is, a linear kernel through 
which stimuli are filtered before they are further processed. If information transmission is assumed 
to proceed in a feedforward cascade, the receptive field may be interpreted as the external stimulus’ 
profile maximizing neuronal output. The nervous system, however, contains many feedback loops, 
and sensory neurons filter more currents than the ones representing the transduced external stim¬ 
ulus. Some of the additional currents are generated by the output activity of the neuron itself, and 
therefore constitute feedback signals. By means of a time-frequency analysis of the input/output 
transformation, here we show how feedback modifies the receptive held. The model is applicable 
to various types of feedback processes, from spike-triggered intrinsic conductances to inhibitory 
synaptic inputs from nearby neurons. We distinguish between the intrinsic receptive held (filtering 
all input currents) and the effective receptive held (hltering only external stimuli). Whereas the 
intrinsic receptive held summarizes the biophysical properties of the neuron associated to sub¬ 
threshold integration and spike generation, only the effective receptive held can be interpreted as 
the external stimulus’ prohle maximizing neuronal output. We demonstrate that spike-triggered 
feedback shifts low-pass hltering towards band-pass processing, transforming integrator neurons 
into resonators. For strong feedback, a sharp resonance in the spectral neuronal selectivity may 
appear. Our results provide a unified framework to interpret a collection of previous experimen¬ 
tal studies where specihc feedback mechanisms were shown to modify the hltering properties of 
neurons. 
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I. INTRODUCTION 


Dong and Atick 1995, 

Schwartz and Simoncelli 2001, 

Geisler 2008, 

Thcunissen et al. 20001 


The dynamic range of input signals can often span several orders of magnitude during 
the course of a single behaviorally relevant time interval. For example, in the visual sys¬ 
tem, mean luminosity changes drastically when the gaze is displaced from a spot that is 
directly illuminated by the sun, to a shadowy corner. Neural systems have therefore de- 
veloped adaptive mechanisms, modifying the neural code according to the sensory context 


Several types of adaptive mechanisms exist, as for example, 

synaptic plasticity 

Atwood and 

Carunanithi 2002, 

Xu-Friedman and Regehr 2004, 

Siostrom et al. 2008. 

Feldman 2009 

, feedback through recurrent connectivity 

Douglas et al. 1995. 

Shu et al. 200, 

3, Eyt 

an et al. 2001 

>. Buonomano and Maass 20091, feedback through 


adaptation currents 


Wang 1998 


, ISanchez-Vives et al. 2000 


Peron and Gabbiani 20091] . and intrinsic non-linear effects 


Prescott and Seinowski 2008 


Borst et ah 2005]. Different 


mechanisms operate on different timescales; whereas non-linear mechanisms emerge in a 
matter of milliseconds, synaptic plasticity usually develops in one or a few seconds. Changes 
in the intermediate range (hundreds of milliseconds) are mostly due to adaptation currents 
and recurrent connectivity. At these time scales, adaptivity arises from the dynamics 
associated to certain processes at cellular and network levels, without the profound 
reorganization entailed by learning and plasticity. 

Adaptive phenomena mediated by intrinsic currents and feedback network connec¬ 


tivity are based on history-dependent spike-evo ked ac tivity. 


back influence mainly by reducing neuronal gain 


Felsen et al. 2002 


l ey exert t heir feed- 


Avaz and Chance 2009i . 


Benda and Herz 20031 . 


Benda et al. 2010], and thereby, by modifying the input/output re¬ 


lation of the cell, that is, the relation between spiking rate and mean stimulus strength. 
However, adaptation phenomena go far beyond a mere reduction in firing rates, often in¬ 
volving a dramatic reshaping of the selectivity to time-dependent inputs, and of the statis¬ 
tics of neuronal output. For example, in the olfactory bulb, feedback has been shown to 
amplify input fluctuations of a specific frequency, and thereby to induce strong oscilla¬ 


tions in the output activity Freeman 1972al. 


Freeman 1972b 


Fr eeman 1 972b ], giving rise to 
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complex (often chaotic) dynamic behavior Freeman 19871. iDavid and Friston 20031 ] . Simi- 


tions 

Wilson and Cowan 1972, 

Coombes and Lai rig 2009, 

Buice et al. 2010. 

Bressloff 2012, 

Arnit and Brunei 1997. 

Ledoux and Brunei 2011 

. Therefore, although adaptation is usually 


claimed to have evolved in order to increase the dynamic range of sensory encoding, its effect 
on the temporal properties of the neural code should not be overlooked: Feedback alters the 
basic properties of neuronal selectivity (as explored below), and also the temporal evolution 
of the out put, the amount of temporal correlations, the precision and th e reliability of neural 


responses 


Ladenbauer et al. 2014, 


Urdapilleta 2011 


Butts et ah 2011]. 


In the context of spike-evoked feedback, adaptation processes have been shown 
to modify the selectivity to transient stim ulus temporal profiles , enhancing the rep- 


resentation of high-free 

mency component 

Benda and Hennig 2008, 

Benda et al. 2010 


Benda and Herz 2003, 


Gigante et ah 200 


ing the recep tive held, or the relevant st imulus directions, by means of reverse correlation 
techniques Dayan and Abbott 200ll . IChichilniskv 200U ISamengo and Gollisch 2013 ]. 
Several previous studies have demonstrated that adaptation in mean hr- 


ing rates may or mav not be accomi 

rallied by changes in receptive helds 

Enroth-Cugell and Shapley 1973, 

Victor 1987 


Butts et al. 2011. 

Baccus and Meister 200^. 

Sharpee et al. 2011. 

Garvert and Gollisch 2013 

j. To our knowledge, there is yet no theo- 


retical framework that allows us to understand why and when receptive fields are expected 
to be modulated by feedback. Such theoretical framework should be general enough 
to be applicable to cases where adaptation is mediated by feedback at the network 
level, or at the level of voltage-dependent_ ionic cond uctanc es and refractoriness, as for 


example, in_the Hodgkin-1 


uxley model 


Samengo et ah 2013] or in a LGN model neuron 


Gaudr v and Reinagel 2007], In this work we describe the adaptive changes observed in 


the stationary and transient encoding properties of linear Poisson models driven by a 
combination of external stimuli and spike-triggered negative feedback. We theoretically 
analyze how adaptation modifies the shape of the receptive held, providing a complete 
spectral and temporal description in the limit where the feedback signal is fully determined 
by the spiking probability of the cell under study (perfect feedback). Spike-evoked negative 
feedback is shown to induce divisive gain control, to reshape the receptive held, and to 
enhance resonant properties. In order to extend these results to the case where feedback 
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is a noisy function of the spiking probability, we incorporate stochastic elements into the 
theory. Finally, we also discuss an extension to non-linear Poisson models. 


II. RESULTS 

A. Theoretical description of perfect feedback 

The main goal of this work is to study how the receptive held of a neuron varies, when 
feedback processes are incorporated. In addition to the temporal dimension, receptive fields 
may be defined as a function of a variety of additional dimensions (spatial, frequency, chro¬ 
matic, chemical) depending on the modality of the sensory system under study (vision, 
audition, taste, etc.). Since adaptation processes unfold along the temporal domain, for the 
moment, we restrict the analysis to the temporal profile of receptive fields, and defer to the 
last section of result the extension to higher-dimensional problems. 

Strictly speaking, the concept of receptive held is well defined for linear-nonlinear Poisson 
models. In order to develop the theoretical framework, we initially restrict to purely linear 
Poisson models, and later on discuss the extension to the nonlinear case. In the linear case, 
the probability of generating a spike in the interval [t, t + dt] under the influence of stimulus 
/(f) is r(t) d t, with 


r(t ) = h 0 + h{r) I(t — r) dr. 


( 1 ) 


Here, h 0 is the spontaneous bring probability, and h(r) is the receptive held of the 
cell. Since only past stimuli can influence the present firing probability, causality im¬ 
poses that h{r ) = 0, for all r < 0. The stimulus I{t) can be interpreted as either the 
external signal controlled by the experimentalist (light, sound, touch, etc.), or as the in¬ 
put ionic current entering the cell. By stimulating the cell with stochastic input sig- 
nals- the shape o f h(r) may be easily obtain ed through reverse correlation techniques 


Chichilnisky 2001 


Sa me ngo and Gollisc h 2 013]. The shape, units and dimensionality of 


the hlter h{r) naturally depend on whether the input signal used in the reverse correla¬ 
tion analysis is the macroscopic external stimulus, or the microscopic ionic current. We 
specifically distinguish between the component of I(t) representing the transformation of 
the external signal accomplished by upstream neurons, and the component describing all 
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negative feedback signals that are triggered by previous activity of the cell under study 


I(t) = s 0 +si(t) - g x(t). (2) 

In Eq. ([2]), So + Si(t) is the extrinsic stimulus component stemming from the transduc¬ 
tion of sensory signals, and x[t ) is the feedback signal whose value depends on the spiking 
history of our neuron. The coupling constant g has dimensions of transduced stimulus, and 
represents the strength of the feedback connection. It may be either negative or positive, 
depending on whether the cell under study fires in response to positive or negative stimulus 
deflections. The sign must be chosen in such a way as to produce a signal that opposes the 
natural excitability of the cell, in order to avoid positive-feedback instabilities (see below). 
We separate the average external stimulus So, so that the time-dependent component Si(t ) 
can be assumed to have zero mean. Defining the baseline firing rate 


r 0 = h 0 + H s 0 , 

and replacing Eq. ([2]) in Eq. ([I]), we obtain 


where 


r(t) 



[si(f t) g x(t 


r)] dr, 


(3) 

(4) 


H= h(r) dr = V2ir h(cu = 0). (5) 

Jo 

The parameter H, hence, is proportional to the zero-frequency component of the Fourier 
transform of the receptive held, h(co = 0) (the convention used for the Fourier transform is 
specified in Appendix 13}. 

We assume that the feedback signal x(t) is boosted by discrete pulses, and has a natural 
decay time iq. For example, if feedback is implemented through the inhibitory action pro¬ 
duced by nearby neurons, the pulses represent spikes produced by other cells in the network 
providing inhibitory feedback to the neuron under study. Hence, 


da: 

df 


N 

ai 

Ta z ' 

u i wi. 


( 6 ) 
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where a x weighs the increase of feedback activity due to a spike in the i-th presynaptic 
neuron, occurred at time tf (k collectively represents all spike times). If Eq. ([6]) is meant to 
represent a feedback signal, the average activity of the neurons contributing to the sum must 
be proportional to the output of the neuron under study. Only with such proportionality can 
we ensure that x(t) is linked to the past activity of the neuron. For the sake of simplicity, we 
model the link as a simple proportionality. Moreover, in our first attempt to model feedback, 
we assume that the sum in Eq. (J6]) is not proportional to the actual output of the cells, but 
rather, to the probability r(t) to generate a given output. That is, 


N 


ai ( 7 ) 

i=l, Vk 

This approximation is here called a perfect feedback signal, since x(t) is a determinis¬ 
tic function of r(t); more specifically, it is a leakily integrated copy of r(t). Later on we 
discuss the case where the feedback signal is a stochastic (as opposed to deterministic) 
function of r(t), and more accurately describes the actual output of the cell under study. 
The approximation of perfect feedback is valid if there is a large number of neurons con¬ 
tributing to the sum in Eq. (J6]) , and if all of them have similar statistical and dynamical 
properties, so that they are all governed by the same firing probability. Alternatively, one 
may assume that the external input drives a pool of independent neurons, which collcc- 


tivelv provides a normalizing signal associated to the common processing 

stimuli 

Heeger 1992, 

Carandini et al. 1997. 

Carandini and Heeger 2012 

]• 


these conditions mean to assume that all the cq are equal, to take N —>■ oo, and in order to 
maintain the total input bounded, additionally scale = 1/N. With these approximations, 
replacing Eq. (JTJ) in Eq. ([6]), we obtain 


dx 

df 


x , , 

-b r(t). 

Td 


( 8 ) 


Equations (J4|) and (JH]) constitute a closed set: If the external stimulus s 0 + Si(t) an d the 
filter h{r ) are known, both the feedback term and the firing probability can be calculated. 

The aim of this study is to deduce how feedback affects the filtering characteristics of 
the neuron. To that end, we now assume that both the external stimulus and the firing 
probability are known (the latter can be recorded from repeated presentations of the same 
stimulus), and evaluate whether the input/output relation can still be considered a filtering 
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FIG. 1: Spectral characteristics of simplified ON, OFF, and biphasic filters. A: Example of an 
ON filter displayed as a function of time (left). The modulus and phase of its Fourier transform 
are displayed as a function of frequency (right). Note the logarithmic and semi-logarithmic scales. 
B, C, D: Analogous representations of OFF, biphasic-ON and biphasic-OFF filters, respectively. 
In these examples, ON and OFF filters have exactly the same \h(uj)\, and the phases are shifted in 
7T (compare A and B ). The same relation is found between biphasic-ON and biphasic-OFF filters 
(compare C and D ). 


process, in spite of feedback. The linear nature of Eqs. ((4j) and ([8]) calls for a treatment 
in Fourier space. Based on the properties of the Fourier transform (see Appendix 0 and 
rearranging terms, we obtain 


r(co) = 


x/27T (1 + i U! Td) 


1 + i u Td + \/2vr g Td h(u) 

x[r 0 S(oj) + h(ui) si (a;)]. 


In the absence of feedback (g = 0), this expression reduces to 


r 0 5(u) + h(u) §i(u) 


(9) 


r(oj) = \/2 tt 


( 10 ) 






























Solving Eq. (ITOj) for the filter at non-zero frequencies leads to 


h(u) = (1/V2 tt) r(oj)/si(oj) 

= (l/v 7 ^) f(cn)s*(a;)/|si(a;)| 2 , 


( 11 ) 


so for white-noise stimuli , the filter in the temporal domain h{r) is proportional to the 
spike-triggered average Gabbiani and Koch 1998]. 

When feedback is active (g ^ 0), the filter h{r ) can no longer be calculated with Eq. (HI]) . 
Comparing Eqs. ([9]) and (fTOlh we see that for each frequency oa, the firing probability is still a 
linear function of the applied stimulus. The constant of proportionality, however, is modified 
due to feedback, and the modification affects differently the continuous component (oj = 0) 
and the non-zero frequencies (uj 0). For the continuous case, feedback changes the baseline 
bring rate ro to an effective value 


fb h 0 + H s 0 _ ^ 


r o = 


= K + H s 0 , 


( 12 ) 


1 + g r d H 

which implements a divisive gain rescaling. 

For non-zero frequencies, feedback modifies the coding properties of the neuron in such 
a way that the intrinsic filter /r(r), processing both the external stimulus and the feedback 
signal, is equivalent to an effective filter hr°{T) that only filters the external signal. In Fourier 
space, the relation between the intrinsic and the effective filters is 


hP 3 ^) =- - + --- h(cu). (13) 

1 + i u Td + \/27r g t d h(u) 

In the remaining part of the present section, we analyze the effects of Eqs. (fT2l) and 
(TT3D both in the frequency and the temporal domains, for several types of filters. We then 
extend the analysis to imperfect feedback processes, where the sum in Eq. (J6J) is no longer 
deterministically proportional to r(t). Finally, we consider non-linear Poisson neurons, and 
we discuss the validity of the linear approximation. 


1. Poisson neuron models in the absence of feedback 

In order to understand the effect of feedback, we first describe the basic types of processing 
in a purely feedforward model. Typically, neuronal filtering characteristics are classified 
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FIG. 2: Temporal evolution of the firing probability of different cells, with and without feedback. 
A-D : ON, OFF, biphasic ON, and biphasic OFF cell models, respectively (see Fig. [I]). Bottom 
panel: Applied stimulus. At time t = 0 ms, a step of amplitude s o = 0.05 stimulus units is 
applied. At time t = 500 ms, the sinusoidal stimulation begins (frequency 10 Hz and amplitude 
0.01 stimulus units). Intermediate panel: Firing probability in the absence of feedback. Top panel: 
Firing probability r(t) (continuous red line; scale on the left margin) and feedback process x(t) 
(continuous gray line; scale on the right margin) for the model with negative feedback. Shared 
parameters: g = 0.005 stimulus units, rj = 100 ms. Monophasic filters are characterized by \H\ = 
2.506 (stim) -1 .(ms)~ , whereas biphasic filters are symmetric, H = 0. In all cases, ho = 300 Hz and 
the response is measured in (ms) -1 . For OFF cells (both monophasic and biphasic) the coupling 
factor g is negative. 


according to the shape of h{r). For example, in visual areas, the filters of simple cells 
are classified in a limited number of types jSegev et ah 20061 ]: ON, OFF, biphasic ON and 
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biphasic OFF cells, as shown in Fig. |U In a previous study Urd apilleta and Samengo 2009], 
we demonstrated that when processing slow stimuli, the firing probability of these four 
different cells is a simple function of the external stimulus. In ON and OFF cells, r(t) is 
proportional to Si(t — 5), where the delay 5 is determined by the shape of h{r). In biphasic 
cells, r(t) is proportional to s', (t — 5), where is the temporal derivative of the stimulus. 
In Appendix |Bj we offer a novel derivation of these results, based on the Fourier approach 
developed in this paper. 


2. Poisson neuron models with perfect feedback 

To study the effects of feedback, we separate the analysis in two: First, we focus on the 
mean response, determined by the spectral content of the response at cu = 0. Later, we 
analyze non-zero frequencies, w ^ 0. 

Effect of feedback on the baseline firing level 

As stated in Eq. (TT2|) . the presence of feedback reduces divisively the firing probability. 
In Fig. [2] we observe how the magnitude of the reduction depends on the type of cell. In 
the time interval between 0 and 500 ms, the response to a constant stimulus is displayed. 
Before the onset of the step stimulus (t < 0), the presence of feedback (g ^ 0, red curve, 
upper panels) reduces the steady state firing probability, as compared to the value obtained 
in the absence of feedback (g = 0, middle panel), for ON cells (Fig. [2]A)- The same holds 
for the steady-state bring probability after the step increase in input current, in the interval 
between ~ 250 and 500 ms. The stationary bring level of OFF cells is also diminished 
(Fig. [2]B) by feedback. In this case, negative feedback is implemented through a negative 
coupling factor, g < 0, which could, for example, correspond to excitatory feedback from 
lateral OFF cells. An OFF cell with a positive g value could cause the bring probability 
to grow unboundedly, and the Fourier transforms to become ill-dehned. For biphasic 
symmetric cells, the integral of the biter is zero, H = 0, so feedback does not modify the 
asymptotic processing of stationary signals (see Figs. [2]C and[2p, for f <0). 

Effect of feedback on temporal processing 

Combining the ebective baseline level of Eq. (fTiZjl and the ebective biter h^iuj) debned 
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in Eq. (fl3|) . we may now re-write Eq. (J9j) as 


r[u 


= \/2tt 


5(o;) + h tb (uj ) Si(o;) 


fb/ 


(14) 


which is formally equal to Eq. (HU . As a consequence, in the time domain, 


r(t ) =r^+ h^{r ) si(t - r) dr. (15) 

J — OO 

The intrinsic filter /i(r) characterizes the biophysical properties of the cell under study. 
Since h(r) hlters both the external component of the stimulus So + si(t) and the feedback 
signal x(t), its shape cannot be calculated through reverse correlation analysis performed 
with solely the external signal; knowledge of the internal signal is also required. The effective 
filter /^(t), instead, can be obtained with the sole knowledge of s\(t). Once h^iuj) is 
obtained, the intrinsic filter h(oj) can be recovered with Eq. (fT3|) . assuming that feedback 
properties are known. The temporal profile of the effective filter h^{r) can be obtained 
simply by inverse Fourier transform. 

The Fourier spectrum of the intrinsic hlters of both monophasic and biphasic cells decays 
at large frequencies (see Fig. [[]). Hence, in the high frequency range, the denominator in 
Eq. (fT3l) is approximately equal to 1 + i co iq, and consequently, lim ^oo |/r b (a;)| ~ \h(co)\. 
Biphasic hlters also have a reduced spectral content in the low frequency range (see Figs. [DC 
and HP). Therefore, in these cells, feedback exerts a limited effect in the whole frequency 
range. Since the filtering properties of these cells are hardly modified by feedback, hereafter 
we focus on monophasic cells (ON and OFF). 

In Fig. HI we discuss the differences in the spectral domain between hP^^uj) and h(oj). For 
any OFF hlter there is another ON hlter that is exactly equal in shape, but with inverted 
sign. To consider negative feedback, the coupling constant g acting on an OFF cell also has 
to be inverted in sign, otherwise, we fall on the case of positive feedback, which is unstable 
in the linear case. Hence, due to inversion symmetry, ON and OFF cells present exactly the 
same gain spectral behavior. Additionally, hlter phases of ON and OFF cells display the 
same spectral characteristics but shifted in n. 

Overall, for monophasic cells, we observe a moderate spectral reshaping of hlters due 
to the presence of feedback. The low-pass filtering characteristic of h(oj) are converted to 
band-pass hltering properties in hP°(oS) (see black and red continuous lines, respectively, in 
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A 




FIG. 3: Spectral composition of filters of monophasic cells, with and without feedback (hP° in red 
and h in black lines, respectively). A: The modulus |/i(u;)| is exactly the same for ON and OFF 
cells. Feedback reduces the spectral content at low frequencies and induces a relative increase in 
the intermediate range. B: Phase of and h(oj), for an ON cell (full lines) and an OFF cell 

(dashed lines). Feedback advances the phases at intermediate frequencies. Parameters as in Fig. [2l 


Fig. EK)- Feedback reduces the contribution of low frequency stimuli and slightly enhances 
the influence of stimuli with intermediate frequencies. Since the system is linear, the ap¬ 
plication of a sinusoidal stimulus (see Fig. [2j for t > 500 ms) evokes a sinusoidal response 
probability; the ratio between relative amplitudes of the response and the stimulus (multi¬ 
plied by a factor l/\/2n) defines the gain or the modulus of the filter. The introduction of 
negative feedback produces a minor effect on the phase. Importantly, however, a phase ad¬ 
vance is observed at intermediate frequencies, implying that r overtakes s in the asymptotic 
regime (see Fig. [2] for t ^ 500 ms). 

Feedback is defined by two parameters: the time constant Td that determines the temporal 
development of x(t), and its relative contribution to the input signal, given by g. In Figs.[4]A 
and |4]C we study the effect of varying the coupling strength g on a monophasic ON filter. 
When g is small, feedback barely influences neural processing. As g increases, low frequency 
content decreases in order to satisfy h fb ( ui) —> i? fb /v / 2vr = h( 0)/(l + g Td H) and band-pass 


13 
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B 






FIG. 4: Influence of parameters defining feedback on spectral and temporal processing. Black 
line: intrinsic filter h. A: Modulus of the feedback filter hr for different values of the coupling 
strength g (red lines), with g = 1, 2, 5, 10 and 20 x 10~ 3 stimulus units. For all cases, ra = 
100 ms. B: Modulus of the feedback filter for different time constants rj (red lines), with 

7d = 20, 40, 100, 200 and 400 ms. For all cases, g = 5 x 10 -3 stimulus units. C, D: Temporal 
filter reconstructed from the spectra in A and B and the corresponding phases (not shown). C: 
As g increases, the region of inverted polarity becomes more prominent. D: As Ta increases, the 
recovery from the inverted polarity becomes slower. In all cases, H = 2.506 (stim) -1 .(ms) . 


behavior is emphasized. To see how these characteristics appear in the temporal domain, we 
calculate the filter hF°(r) by applying the inverse Fourier transform on hF°{oS). The effective 
filter /r fb (r) is no longer purely monophasic, as the original filter h(r), since it contains a late 
phase of inverted polarity (see Fig. WP)- The absolute value of the integral of the filter is also 
reduced, when compared to the original H. Therefore, the feedback filter is more sensitiv e 
to the fluctuations of the stimulus than the original filter Urda pilleta and Samen g o 200 9], 
As the coupling strength g grows, stronger feedback generates a more significant region of 
inverted polarity in /i fb (r). 

The effect of varying Ta is shown in Figs. H]B and 0p. For fixed g, if increases, low 
frequency components decrease (Fig. S]B). In addition, for small values of r c i, the gain extends 
its zero-frequency value, H®°/ \/27r, into a wider range of positive frequencies. As shown in 
Fig. [4]B, Ta has only a minor influence at intermediate frequencies. Finally, in the temporal 
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domain, the duration of the region of inverted polarity is governed by r d (see Fig. [4p). 

To summarize, we conclude that whenever reverse correlation reveals a filter with biphasic 
characteristics, the obtained effective filter /^(r) may not coincide with the intrinsic filter 
h{r). In particular, at least some of the biphasic filtering characteristics may derive from 
negative feedback. 

3. Feedback-induced resonances 

A monophasic receptive field behaves as a low-pass filter, with a cutoff frequency deter¬ 
mined by the inverse of the duration of the non-zero portion of h(r). The effect of feedback 
is to additionally reduce the response at low frequencies, so the effective filter h fb (a;) ac¬ 
quires band-pass characteristics. The value of the lower cutoff frequency and the resulting 
Q-factor depend on the properties of feedback (g and r d ). It is therefore conceivable that by 
appropriately choosing these two parameters, feedback can be shaped as to induce a strong 
resonance in the system, even up to the point of instability. In Eq. (I13p . this kind of strong 
resonance appears as a sharp peak in h fb (o;), or even a divergence. The denominator of 
Eq. (1I3D can indeed vanish for specific combinations of discrete frequencies ujj and values of 
the product (g r d )j. In these cases, the presence of negative feedback renders the system un¬ 
stable, amplifying one particular frequency (or several). Near a resonance, any infinitesimal 
stimulus component matching the critical frequency is amplihed in the response, giving rise 
to strong oscillations. The oscillations enter repeatedly into the feedback loop producing a 
divergent response. Of course, no real neuron can truly produce diverging responses, because 
as oscillations grow in amplitude, Eq. (0Q) loses validity: The evolution of the system can no 
longer be described by a linear equation; in particular, the linear scheme must be abandoned 
before the oscillations in the firing probability r{t ) are strong enough as to produce negative 
values. 

The linear analysis is nevertheless useful to point out the dramatic amplifying effect that 
negative feedback can have, and the conditions that favor resonances. Clearly, if the denom¬ 
inator of Eq. (U3D vanishes for one or more frequencies, the resulting effective filter h^ioj) 
cannot be transformed back into the temporal domain. As stated above, such divergences 
may appear (if at all) at a discrete collection of frequencies ujj, and discrete values of the 
product (g r d )j. The absence of divergences, however, does not guarantee that the inverse 
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FIG. 5: Feedback-induced resonances for the ON filter of Fig. |T| A: Modulus of the effective 
filter |/r fb (w)| near the resonant instability. For negative feedback with timescale ra = 100 ms, a 
sharp peak develops at w « 320.4 rad/s. In this case, the critical value for the coupling strength is 
gi ~ 0.1346 stimulus units. B: The phase of the effective filter varies rapidly near the resonance. 
C: In the temporal domain, as the resonance is approached, the effective filter develops strong 
oscillations at the frequency of the critical cjj. For g > gi, the filter h™(u) cannot be transformed 
back to the time domain, implying that the linear system of Eq. m is ill defined. 


Fourier transforms hF a {r) and r(t) exist and are bounded. In fact, the condition that gives 
rise to instabilities is broader and includes the discrete cases where the denominator of 
Eq. (1T3l) vanishes. The proper ma thematical framework to analyze the onset of instabili¬ 


ties is provided by control theory Franklin et ah 1994 ], and can be addressed in terms of 
the behavior of the Laplace transform of the effective filter h fb (s) with complex argument 
s — a + i uj. As derived in Appendix [A] an instability appears when at least one pole 
s j = CTj + i ujj has positive a,,. The pole s* with largest real part is hence critical, since the 
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FIG. 6: Inhibitory feedback supported by a finite population of neurons. The signal s\(t) stimu¬ 
lates a population of identical neurons. The output spikes are filtered during synaptic transmission 
and then return in the form of a negative feedback input current. 

magnitude of cr* determines the stability of the system. When a* < 0, the effective filter 
is qualitatively similar to the examples shown in Figs. OandHJ As the coupling strength g 
increases, s* gradually shifts to the right and gets closer to the imaginary axis; consequently, 
the effective filter begins to show a prominent peak at the frequency oj* of the critical pole, 
as shown in Fig. [5]A- This resonance is caused by feedback. 

At the resonance, to ~ c o*, the phase of /r fb (o;) varies rapidly (Fig. [5]B). In the temporal 
domain, the effective filter exhibits strong oscillations (Fig. EJC), which occasionally grow as 
far as to make the system unstable. A mathematical analysis of the conditions giving rise to 
instability (see Appendix O reveals that unstable behavior is only observed if g is above a 
critical threshold, the value of which is determined by the shape of the intrinsic filter h(oj). 
For large g, however, the baseline firing level r q drops significantly (see Eq. (fT2jl ). Therefore, 
a transition to instability may only be expected in systems with strong feedback and, si¬ 
multaneously, with large intrinsic spontaneous activity ro or decreased stimulus fluctuations 
,s'i (f), so as to ensure that the firing rate r(t) remains positive. Once these conditions are 
met, the frequency to* of the unstable oscillations is determined by the time constant iq and 
the shape of the intrinsic filter h[r). In the limit of large ra (as in the example of Fig. [5]), 
the feedback time constant becomes irrelevant, and the location of the critical frequency 
depends only on h(r). 
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B. Poisson neuron models with imperfect feedback 


So far, the analysis was based on Eqs. ([4j) and (ED, where the temporal derivative of 
the feedback signal x(t ) was proportional to the bring probability r(f). This relation is only 
valid in the limit of a homogeneous population of infinitely many identical neurons, uniformly 
coupled and driven with identical stimuli s(t). These conditions are hardly realistic, since 
the amount of feedback must be determined by the actual number of generated spikes, and 
not by the spiking probability r(t). There is an important difference between these two 
options. The spiking probability r(t ) is a deterministic function of the input current J(f) 
(Eq. ©), whereas actual spikes are stochastic point processes governed by r(t). Therefore, 
in an attempt to provide a more realistic description of feedback, we now assume that x(t ) is 
given by Eq. (ED . Although this new model can be studied analytically in specific parameter 
regimes (see below), we initially resort to a numerical approach to determine the point up 
to which the results of the previous section can be extended to more realistic conditions. 

In the present description, illustrated in Fig. El we assume that neurons are limited in 
number, all have identical intrinsic filtering properties h(r ), and process a common temporal 
stimulus, Si(t). Furthermore, the feedback signal x{t) is assumed to be the same for all cells, 
so the bring probability r(f) given by Eq. (JED holds for any cell in the population. However, 
since spike generation in Poisson processes is stochastic, the precise temporal location of 
spikes dibers from neuron to neuron and, therefore, the sum in Eq. (ED is a random variable, 
which only recovers its deterministic limit (Eq. (JED) when the number of neurons tends to 
inbnity [N —> oo). Since the derivative of x(t) is no longer strictly proportional to r(f), 
feedback is now called imperfect. 

In Fig. [7]A we show the evolution of the bring probability r(t) (upper panel), along with 
the feedback signal x{t) (middle panel) and the driving stimulus (bottom panel) for the 
bnite population model. These three signals are common to all neurons. The feedback 
signal is the sum of N bltered spike trains, a few of which are displayed in Fig. [7]E>. The 
discontinuities in the traces (inset) are produced by individual spikes. These irregular traces 
barely resemble the deterministic counterpart XdetH)/-^ obtained in the limit N —> oo (black 
curve). Their sum, however, smoothes huctuations out, and follows Xdetft) closely (middle 
panel of Fig. EH)- 

To quantify the stochastic behavior of this system, we drive the cells with a sinusoidal 
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FIG. 7: Evolution of dynamic variables in the presence of imperfect feedback. A: Time-dependent 
firing probability r(t ) and feedback process x(t ) for the model with imperfect feedback driven 
by the stimulus s o + si(t). Black lines: prediction for N —>• oo. Irregular red lines: single-neuron 
signals, obtained from a population of N = 10 neurons. B: Single spike trains filtered with a simple 
synaptic dynamics, dx;/d t = — x\/t .a + (1/N)5(t — tf ). Individual signals in colored irregular lines. 
Irregularity arises from the discontinuities evoked by incoming spikes, amplified in the inset. The 
feedback signal x(t) shown in A (middle panel) is YliLi x i(t) (notice the scale difference with Xi(t) 
in B). Black line: proportional contribution of a single X{(t) in the limit N —>• oo. C: Extracting the 
relation between stimulus and response in noisy conditions. For a linear system, a small sinusoidal 
stimulus (black line; scale of the left) produces a sinusoidal firing probability of identical frequency 
(gray line, scale on the right). A long run of the output spike train is split in T-windows, where 
T is the period of the sinusoidal input. All spikes are then wrapped in a single window, and the 
histogram r(t) is constructed. Top/bottom panels: Stimulation at 10 Hz / 100 Hz. D: From the 
fit of the histograms in C, we calculate the gain and the phase, <j) r — (j) s , of the transfer 

function (Ar is the amplitude of the sinusoidal function that best fits the response). Different 
symbols represent different population sizes. Spikes are collected in asymptotic conditions (that 
is, after the initial transient) during 500 s. For each data point, 10 repetitions are simulated. 
Error bars indicate the standard deviation of the gain and phase, as calculated from the best 
fit parameters obtained for each of these repetitions. Parameters: H = 2.506 (stim) _1 .(ms) , 

Td = 100 ms, g = 0.005 stimulus units, ro = 0.3 mT 1 , sq = 0.05 and As = 0.005 stimulus units. 
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stimulus. As observed in Fig. [7]A, for times t > 500 ms, the firing probability has a marked 
periodic component at the frequency of the input signal. Therefore, the input /output proper¬ 
ties can be characterized by studying a single stimulus cycle. This reduction is implemented 
by taking the spikes fired in different time windows (one window per stimulus period T) 
and wrapping them together into a single window, taking care of preserving the original 
firing phase with respect to the stimulus. That is, each spike is displaced an integer num¬ 
ber of periods, and located within a window whose duration is equal to a single period of 
stimulation. Once a long spike train realization is so wrapped, the corresponding histogram 
in the T -window is built, as shown in Fig. [T]C. In this figure, we can clearly observe the 
periodic modulation of the response (stairs-like gray lines), which expectedly can be fitted 
by a sinusoidal function (continuous blue lines). The amplitude and phase (relative to the 
stimulus) of the adjusted response is used to construct the spectral characteristics of the 
stochastic model. 

In Fig. 0C we show the spectral characteristics of the population-based feedback model. 
Both the gain and the phase shift are practically invariant with the number of neurons in 
the population, and a good agreement with the case N —* oo (black line) is observed. The 
dispersion at each point (error bars) arises from the variability in the sinusoidal fit, which 
depends on the irregularities of the histogram (controlled by the total recording time and the 
mean number of spikes produced during a cycle). The irregularities, in turn, arise from the 
Poissonian character of the spike generation process. When N = 1, the filtered activity of the 
neuron itself is used as the feedback signal defined in Eq. (J6J), and the sum in % just involves 
a single element, i — 1. This situation is adequate to model self-inhib i tion due to spike - 


triggered adaptation currents 


Benda and Herz 2003 


Benda et al. 2010, U rda pilleta 2011], 


As N grows, the description gradually shifts to represent network-mediated feedback pro¬ 
cesses. 

The wrapping procedure used to construct Fig. |T| allows us to determine the linear re¬ 
sponse function as the ratio between the amplitude of the firing probability (response) and 
the stimulus driving the system (input), both measured at the same frequency. One impor¬ 
tant result is that the linear response function is independent of the population size N (see 
Fig. EP). However, the irregularities observed in the temporal response do indeed depend 
on N. To exemplify this behavior, in Figs. [8]A and [8]B we show the time-dependent firing 
probability r(t ) and the feedback process x(t), for the cases N = 10, N = 5 and N = 1. 
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Irregularities become more prominent as the population size decreases. When the number 
of neurons is finite, the firing probability contains a certain amount of power at frequencies 
different from the incident frequency. The power at spurious frequencies does not affect 
the linear response function, because the fitting procedure is casted specifically at the input 
frequency. So far we have described the properties of the response at this input frequency. 
Now we turn our attention to the rest of the spectrum. 

As an example, in Figs. HP and [8]E we show the power spectral density of r(t) and of 
x(t), when the system is stimulated with a 10 Hz, small-amplitude sinusoidal signal. Given 
the linearity of the system, the firing probability and the feedback signal contain a strong 
component at precisely 10 Hz (see the peaks in the corresponding spectra). The height of 
the peak in r(t) is given by the linear response function studied before (see Fig. [7]). The 
remaining spectral power (the background) arises from the inherent randomness of Poisson 
processes. As the population size increases, fluctuations in the spontaneous regime diminish 
and, correspondingly, the background spectral density decreases as well. 

The remaining part of this section is devoted to obtain an analytical expression of the 
background spectrum. We first analyze the feedback process corresponding to N — 1 (adap¬ 
tation current), and later generalize the result to arbitrary N. For N — 1, the signal x(t) 
evolves according to 


dx 

dt 


x s 

-i- £(£)> 

Td 


(16) 


where £(t) represents the spike train produced by a Poisson process (a realization), with a 
time-dependent rate given by Eq. 03]). This rate is coupled to Eq. (1TB]) (and thereby, to the 
noise source) through the feedback term. By approximating the time-dependent bring rate 
by its baseline level r^ (Eq. f|T2|) ). Eq. (fl6|) becomes a standard linearly bltered Poisson 
process (bltered shot noise). Under this approximation, it is simple to hnd the formal 
solution to Eq. f|T6|) and the resulting exponential autocorrelation function, 


1 


Cx(r) = r d e 


-MAd 


(17) 


The power spectral density (one-sided, per unit time) of the feedback signal x(t) 
can be obtained from the autoco rrelation function through the Wiener-Khinchin theorem 


Gardiner 1985 


Press et ah 2007 ], 
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FIG. 8: Characterization of finite size effects. A, B: Evolution of the time-dependent firing rate 
r(t) (in ms^ 1 ) and feedback process x(t) (adimensional), for the stimulus in C. In A and B, different 
panels correspond to different population sizes: N = 10 (blue, top panel), N = 5 (green, middle 
panel), and N = 1 (red, bottom panel). Black lines: deterministic limit N —>• oo. All parameters 
are the same as those used in Fig. [2]A. As N decreases, fluctuations are amplified. D, E: Power 
spectral density per unit time of r(t) ( D ) and x(t) ( E ), in the asymptotic regime of small-amplitude, 
periodic stimulation (frequency: 10 Hz, amplitude: 0.005 stimulus units). Colors as in A and B. 
The theoretically derived power spectral density, under the weak coupling assumption, is shown in 
dashed lines (with corresponding colors) for S r (f ) (D), Eq. (1251) . and S x (f) (E), Eq. (|22l) . F, G: 
Analogous to panels D, E for a weaker feedback, g = 0.001 stimulus units. 


S x (u) 
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[|£(-u;)r + |£M| 2 ] 



C x (—u>) + C x (u) 


Since the Fourier transform of Eq. (TT7T) reads 


(18) 
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C x (u) = C x {-u) = 


„fb 


y/2ir (l/r d ) 2 + tu 2 ’ 
the one-sided power spectral density per unit time is 


S x (u) 

T 


„fb 


7T (1/rd) 2 + w 2 ’ 


( 20 ) 


or, in terms of frequency, 


S x (f) = 2?r S x (u) 


9 ~fb 
Z / Q 


( 21 ) 


T - T (1/^d) 2 + (2?r /) 2 
The previous analysis can be easily extended to the case N > 1. In this case, the effective 

rate for £(t) in Eq. (IT51) is now N r q \ We recall that, in order to maintain the feedback 
level constant when N increases, the efficacy of each spike in Eq. flBJ) has to be proportional 
to 1/N. Therefore, now dx/dt = —x/r d + (1/iV) £(t), and the autocorrelation function is a 


scaled version of Eq. (fT71) , C, 


N neurons/ 


z V r ) = C x 


1 neuron/ 


t)/N. Consequently, 


S*{f) 


9 ».fb 

z r 0 


( 22 ) 


T N (l/r d ) 2 + (2ir/) 2 
This expression is represented with dashed lines in Fig. [8]E, for different population sizes. 

The power spectrum obtained from simulations agrees with the theoretical description, ex¬ 
cept at low frequencies, where higher order statistical interactions between x(t) and r(t) 
become noticeable. At medium and high frequencies, the theoretical approach provides a 
very good description of the simulations. 

The bring probability r(t) inherits the correlation structure of x(t), 


C r {r) =g 2 h(r') dr' / h(r") C^\r + r' — r") dr", 

Jo Jo 


1 

2 N 


g 2 r o T d / h(r') dr' 


h{r" 


X e -\r+r'-r"\/r d 


(23) 


The hlter h[r) is different from zero inside a hnite window [0, r m ] 
Urdapilleta and S amen go 2009]. Therefore, in Eq. (125]) . we replaced the upper inte¬ 


gration limits by r m . Whenever r m -C r d , the exponential factor in the integrand of Eq. 
can be further simplified, so that only time differences r comparable or larger than r, 
matter: 
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• r ~ 0{r m ) => e -l r+T '- r "l/ rd w 1, 


• r» 0(r m ) =>• e l T+r, r "l/ Td « e l r l/ Td 


In this case, the autocorrelation function can be approximated by 



( 24 ) 


Finally, based on this equation and the Wiener-Khinchin theorem, the power spectral 
density for the firing probability r(t) is 



(25) 


This expression is represented in Fig. |8p for different population sizes, in dashed lines. 


As in the previous analysis, at low frequencies, the coupling between x(t) and r(t) produces 
a small discrepancy between the numerical results and the theoretical expression. At high 
frequencies, a faster decay than the predicted S r (f) ~ 1 /p is observed, originated by the 


approximation at r ~ 0(r m ) during the assessment of the autocorrelation function C r (r). 
The derivation in this section is based on the hypothesis that the firing rate r[t) could be 


approximated by its baseline level Tq. This assumption is valid if S\(t) is small, and if x(t) 


and r(t) are weakly coupled. Obviously, as the coupling strength g becomes smaller, the 
theoretical expressions become more accurate, and are valid in a wider range of frequencies. 
As shown in Figs. [S]F and[S]G, in this limit, the simulated spectral densities are in excellent 
agreement with Eqs. (122]) and (123]) . 

C. Linear-nonlinear Poisson neuron models 


Linear Poisson models are only an approximate description of the processes governing 
neuronal dynamics. We here improve the appr oximation by adding a static nonlinearity, as 


often done in the description of sensory systems 

- -7 -— — ' 

Davan and Abbott 2001, Chichilniskv 2001, 

Baccus and Meist< 

3r 2002, Gaudrv and Reinagcl 2007 


Sharnee et al. 2008, 

Butts et al. 2011. 

Garvert and Gollisch 2013]. 

The new model constitutes the linear- 


nonlinear Poisson approach, where 



h(r) [si(t - r) - g x(t - r)] dr 


(26) 


24 
























In this section we demonstrate that although the results of the previous section do not 
strictly hold in the presence of a nonlinearity, it is possible to develop an approximate version 
of the theory that makes very good predictions in most practical cases. 

In Fig. [U]A, we show three sigmoidal nonlinearities /, all with the same functional shape, 
but with different scaling parameters in the x-axis. Here we show how to produce an 
approximate linearized model, and discuss its range of validity. The first step is to find the 
operation point, that is, the spontaneous firing rate Tq 1 , that may differ from the spontaneous 
rate of the linear case r^ (Eq. (TT2]) ). In the absence of time-dependent external stimuli, 
Si(t ) = 0, the linear-nonlinear model reduces to 

= / fa - 9 T d H r ;>] =c /-y; 1 ) =r„-gT d H rg (27) 

Equation (1T7|) implicitly defines Tg 1 as the intersection between the function / -1 (r q 1 ) and 
the straight line r'o — g Td H Tq 1 . In Fig. [9]B, we illustrate the procedure. The level ro (arrow 
on the left margin) sets the offset of the straight line (thick black line). The slope of the line 
is determined by the intrinsic properties of the neuron ( H ) and feedback (g Td). Depending 
on the steepness of the nonlinearity, one same r 0 may elicit different spontaneous bring rates 
rg 1 (colored circles located at the intersections). 

If the constant stimulus component So is modified, a new value ro is established and, 
after a brief transient evolution, the operation level Tq 1 sets to a new value, as predicted 
by Eq. f[2T|) . The y-intercept of the straight line shifts, thus displacing the stationary firing 
probability. In Fig. |9]B the shift is represented by the two parallel lines (thin gray lines), cor¬ 
responding to the addition of two constant stimuli of opposite signs. Due to the non-linear 
nature of the model, the same stimuli displace the response by different amounts (lengths 
of the error bars). Moreover, positive and negative stimuli produce effects of different mag¬ 
nitude (compare the lengths of the left and the right portions of the error bars). 

In Fig. [9]C the temporal evolution of the bring rate (top panel) and the feedback process 
(middle panel) is shown, for a temporally complex input signal. In the spontaneous regime, 
t < 0 ms, the bring probabilities obtained in Fig. [U]B are indicated with colored circles. 
Following the application of a step stimulus, the bring probability and the feedback process 
undergo a transient evolution that rapidly settles onto a new stationary value. When the 
curvature of the nonlinearity is mild, the evolution is similar to the linear case (compare the 
evolution of the red line in Fig. [9]C with the one of Fig. [2]A) . As the nonlinearity becomes 
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FIG. 9: Effects of a static nonlinearity on the feedback model. A: Sigmoidal nonlinearities modeled 
as error functions , (1/2) r max {erf[(g — q c )/ A] + 1}, with A = 0.05, A = 0.10, and A = 0.20 (blue, 
green, and red lines, respectively). In all cases, q c = 0.25 and r max = 0.50. Axis units in ms -1 . B: 
The spontaneous firing rate is defined by the intersection between the line (ro— g t& H r) (thick black 
line) and the inverse function of the nonlinearity (colored curves in A). Spontaneous bring rates for 
the different nonlinearities are indicated by colored circles, for ho = 0.30 ms -1 . The addition of a 
static positive or negative stimulus (gray lines) shifts the operation point nonlinearly (circles plus 
error bars in the upper part of the figure). C: Temporal evolution of the firing probability (top 
panel) and the feedback process (middle panel), under the effect of a particular time-dependent 
stimulus (bottom panel). In the absence of stimulus, t < 0 ms, spontaneous firing rates are defined 
by the colored circles in B. The step stimulus deflects the firing probability transiently, until a new 
stationary value is reached. A periodic stimulus, t > 500 ms, evokes a periodic response probability, 
whose amplitude is a non-linear function of the stimulus’ amplitude. D: Spectral composition of 
filters’ gain in the non-linear model. Continuous black line: gain of the feedback filter for the 
linear model, see Fig. [3]A. Colors match the non-linearities shown in A. Symbols: ratio between the 
amplitudes of the sinusoidal response and stimulus, for a small sinusoidal stimulation, as a function 
of the driving frequency, (1 / y/2n) Ar / As. Irregular lines: gains of the Fourier transform of the 
filters obtained from the spike-triggered average of the nonlinear model, equivalent to x°{u)) for 
small stimuli, see Eq. (|35l) . Dashed lines: equivalent linear filter of the non-linear model /t nl,fb (cu), 
which can be approximated by an appropriate scaling of the linear filter, see Eq. (1371) . Parameters 
are identical to those used before for the linear model (see Figs. [2] and [3]). 
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steeper (blue line in Fig. [9]C), the firing probability evolves faster, and the initial transient 
becomes stronger. 

The next step is to determine how time-dependent stimuli are processed. When a sinu¬ 
soidal signal is applied (t > 500 ms, in Fig. EF), the response probability is periodic, but 
not necessarily sinusoidal. Only in the limit of small s(t) is the sinusoidal response guaran¬ 
teed, and in this limit, the ratio between the (mean-subtracted) relative amplitudes of the 
input/output sinusoidal modulations is independent of the amplitude of the stimulus. The 
amplitude of the response, however, depends on the steepness of the nonlinearity (Fig. EP)- 

To understand the spectral processing of the non-linear model, we define the transformed 
bring rate 


so that Eq. (1261) becomes 


q(t) = f '[>•(«)], 


(28) 


q(t) = r„+ h(r ) [si(( - r) - g x(t - T )] dr. 


(29) 


In terms of q(t), the feedback current x(t), previously described by Eq. (J8J) , is now gov¬ 
erned by 


dx x „ r . . . 

-77 =-+ / [?(*)] • (30) 

dr r d 

Equations (I29T) and (l30j) constitute a closed system, but the presence of the nonlinearity 
precludes the application of the linear Fourier approach that allowed us, in the previous 
section, to find the relation between h{r ) and h fb (r). Linearizing Eq. around the 

operation point, qo = / _1 (rg 1 ), yields 

+ 7 o 1 + f , ( qo )[q(t)-qo]. (31) 

dr r d 

Equations (l29l) and (j3Tj) are closed and linear, and therefore, allow for a linear treatment. 
Transforming them both to Fourier space, and after some algebraic manipulations, 
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q(uj) = y/hr 


q 0 5(uj) 


(1+JW Td) h{oj) 


■silo;) 


Equation 


1 + i uj r d + V2n g r d /'(g 0 ) h(u) 
can be written in terms of an effective filter x fb (o;), such that 

q(u>) = \f2n [®(5(w) + ^V) »i(“)] . 


(32) 


(33) 


where 


x tb M = 


(1 + i u> Td) h{oj) 


1 + i u T d + y/2n g r d /'(g 0 ) h(u) 
The possibility of summarizing the effect of feedback in Eq. 
temporal domain, 


(34) 


implies that in the 


r(t) = f 


"+OO 


Qo + 


X fb ( r ) Si (t - t) dr 


(35) 


The presence of the nonlinearity implies that the effective filter can no longer 

be calculated as the ratio between f(u) and si(o;). Reverse correlation, however, still 
allow us to calculate \ P o {t) from the sp ike triggered average of the recorded data 


Gabbiani and Koch 1998 


Chichilniskv 20.01J. The result is illustrated by the irregular lines 


in Fig. [9p. 

Reverse correlation provides the best possible estimate of / \ fb (cu), since the approach only 
entails the linearization of Eq. (jHJ). The method, however, requires large amounts of data to 
converge to a reliable estimation (2.5 x 10 6 spikes, in Fig. [9p). For practical purposes, hence, 
one may be willing to sacrifice some modeling accuracy, and further linearize Eq. (l35]h for 
the sake of obtaining an easier estimation method. Such approximation brings the response 
and the stimulus to be linearly related 


"+CXD 


r(t) ~ /(go) + /'(go) / sx(t - r) dr, 

J —oo 

so that the linear response function of the non-linear system is 


(36) 











(37) 


_(1 + iu hi) /'(go)_ 

1 + i u r d + y/2n g r d /'(g 0 ) h{u) 
differing from the linear approach (Eq. (If 3D ) by a rescaling controlled by the gain at the 
operation point. 

In Fig. [UP, the spectral behavior of the filter obtained from the spike-triggered average of 
the non-linear model, x fb (a;), as well as the linearized filter h nl,fb (uj) are shown. The symbols 
represent the results obtained from the ratio between f(cu) and si(cj), for small sinusoidal 
stimuli (defined as in Fig[7p). As observed, both filters are in excellent agreement with the 
simulated data. Clearly, as the non-linearity becomes steeper, the band-pass characteristics 
are more pronounced (see also, temporal responses in Fig. [Up). 

The processing differences so far described for nonlinearities of varying steepness are also 
observed in a single nonlinearity, at varying operation points. In Fig. [TOlA. the operation 
point is varied by manipulating the value of h d . The stationary firing rates are represented 
in Fig. [TOl B by symbols located on the nonlinearity (thick black line), whereas the gains are 
indicated on the derivative (gray line). The comparison between the nonlinear model with 
the linear approximations is similar to the comparison of Fig. [Up. 

In the purely linear case, we observed that in certain conditions, feedback could give 
rise to unstable dynamics. Resonances, understood as peaks in the power spectrum of the 
effective filter, are still possible in the nonlinear case. However, now the resonant peak 
cannot grow indefinitely, so unstable divergences are ruled out. The linear approximation 
of Eq. (13TT) has a limited range of validity. For sigmoidal nonlinearities, as the firing rate 
increases, the operation point rg 1 shifts upwards, so the slope /'(go) should not be taken as 
fixed. The value of the derivative diminishes progressively, as the flat part of the nonlinearity 
is approached. In Eq. (134|) . a diminished /'(go) is equivalent to a smaller feedback coupling 
constant g, thereby precluding divergences. Hence, although the firing probability can still 
contain a strong oscillatory component in the nonlinear case, the increased bring at the peaks 
of the oscillations acts as a self-regulatory mechanism, that forestalls unstable dynamics. 

In the present section, we have introduced a nonlinearity in the relationship between r(t) 
and Si(t). The feedback process, however, always remained linear. One could then wonder 
what is the effect of maintaining a linear relation between r(t) and Si(t), but introducing a 
nonlinearity in the feedback process. From Eqs. ([29]) and ([30]) . we see that in terms of g(t), 
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FIG. 10: Influence of the operation point in the nonlinear feedback model. A: A sigmoidal 
nonlinearity operated at different points. Different straight lines correspond to different values 
of ho (0.1, 0.2, 0.3, 0.4 y 0.5 ms” 1 , respectively). The y-intercept is ho + H so, with H = 
2.506 (stim) -1 .(ms) —1 and so = 0.05 stimulus units, and the nonlinearity is defined as in Fig. [9]A, 
with A = 0.10. B: Firing probabilities, for the different values of ho defined in A, indicated by 
symbols located on the nonlinearity (thick black line, scale on the left margin). Associated gains 
are represented with corresponding symbols on the derivative of the function (gray line, scale on the 
right margin). C: Ratio between relative amplitudes of response and stimulus signals, for different 
driving frequencies. Symbols and colors correspond to those represented in A. Black thick line 
displays the spectral behavior of the lineal model, whereas colored dashed lines are appropriately 
scaled versions. 


a nonlinear relation between r(t) and s\(t) traduces into a nonlinear system of equations, 
where the nonlinearity appears in the feedback equation. Conversely, a nonlinearity in the 
feedback process that modifies the second term of the right-hand side of Eq. (JHJ) can be 
treated, from the mathematical point of view, in the same way as the the nonlinear system 
explored in this section. More general nonlinear feedback processes, for example, modifying 
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the first term of the right-hand side of Eq. (J8J), or even mixing the two terms together, 
require additional techniques. 

D. Higher-dimensional receptive fields 

The results obtained so far remain unchanged if the stimulus (and therefore also the 
linear filter) depends on several dimensions, and not just time, as long as adaptation is 
homogeneous in those dimensions. In this case, we write the external stimulus as s(y, t), 
where y G T> is a vector in a one- or a multi-dimensional space T> representing the relevant 
features of the sensory modality under study (direction of the incoming light or chromatic 
composition in the case of vision, frequency content in the case of audition, etc). The 
external stimulus is filtered by some sort of transduction process and perhaps also by one 
or more neurons that lie between the sensory receptors and the neuron under study. We 
represent such upstream filtering processes by means of a high-dimensional filter h u (y,r), 
that transforms the high-dimensional stimulus s(y, t) in a purely temporal ionic input current 


si(t) 


+OO 


K(y, t) s(y,t 


r) dy dr. 


(38) 


J T> J — oo 

This expression can be inserted into Eq. (J2J) to construct the new full input current I(t) 
needed in Eq. (jTJ) . Notice that in the present framework, the high-dimensional stimulus is 
first filtered in the y-dimensional and temporal domains by h u (y,r) (Eq. f[38]l ). and after¬ 
wards by the temporal filter h(r) (Eq. (jTl) ). Using a Fourier analysis completely analogous 
to the one developed before, we arrive at 


r{yj) = 


\f2n (1 + i U! Td) 


1 + i u r d + y/2n g r d h(u) 


r o 


V2n h(u) / h u (y,uj) s(y,u) dy 


IV 


(39) 


This expression reduces to Eq. (J9]) when h u ( y, r) = S(y — yo) <5( r )- 

By stimulating the cell with signals that are localized in the additional dimensions, 
s(y,t) = 5(y — yo) f(t), the whole of the previous theory becomes valid for each cho¬ 
sen y 0 . In particular, feedback still transforms the intrinsic receptive field in an effective 
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receptive field by multiplying the Fourier transform of the former by a factor that depends 
on the frequency, but does not depend on the stimulation point yo- At least, such is the 
effect of feedback if we may assume that the signal x(t) only depends on the output of the 
cell and is not modulated, for example, by spatial input components. 


III. DISCUSSION AND CONCLUSIONS 


In this paper, we analyzed how spike-evoked negative feedback modifies the effective 
receptive held of a cell. The approach was based on an ideal concept, here named perfect 
feedback, where the signal x(t) is not a function of the actual spikes generated by the 
neuron, but rather of the probability that spikes be generated. This assumption is ultimately 
unrealistic, but becomes a good approximation of the real system when (a) the stimulus 
varies slowly compared to the inter-spike interval of the neuron under study, or (b) feedback 
is mediated by a large number of similar neurons in the network. In both cases, the bring 
probability r{t) is sampled exhaustively, and therefore, the distribution of sampled signals 
follows the probability r(t) closely during the time scale Td governing the feedback process. 
This idealized scenario allowed us to develop an analytical approach, and to derive the 
mathematical connection between the intrinsic filter h(r) and the effective filter /i fb (r). In 
particular, Eq. (fT3j) provides the link by which intrinsic or feedback parameters shape the 
receptive held. 

Previous studies have reported some reshaping of receptive helds as the input changes. 
For example, in the case of visual stimuli, the spatial and the temporal context alter the in¬ 
put/output transformation, as the s ystem adapts to the local statistics Schwartz et al. 2007I. 


Wark et al. 200 


1 


War 


Schwartz et al. 2007, 


Kohn 200 


I 


et al. 2009| . producing changes al l the way up to the perceptual level 


Lochmann et al. 2012]. In single cells, both the total lu¬ 


minance and the amount of cont rast s hape linear hlters with incr eased band-pass charac¬ 
teristics 


Enroth-Cugell and Shaplev 1973 


Kaplan and Benardete 2001]. The same effects 


are observed in single audi tory ce lls, as the mean and the variance of synthetic sounds are 


manipulated 


Nagel and Doupe 2006]. According to our results, these phenomena could be 


explained by modulating the amount of feedback as the signal varies. The modulation could 
be mediated by synaptic scaling, in the case of network-based feedback processes, or by ionic 
mechanisms, in the case of single cell adaptation, eg. by Ca 2+ concentration. Moreover, the 
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reshaping of filters can also be modulated by the spatial structure of the stimulus. For 
example, the receptive holds processing stimuli with natural spatio-temporal statistics are 


different from those obtained from simple ensembles. When t 
tics are taken into account in the estimation of receptive fields 


le higher-order input statis- 


Theunissen et ah 2001] . the 


temporal profile of the filters processing natural stimuli turn out to be 


monophasic filters are obtained in responses to gratings David et al. 2004]. In terms of our 


diphasic, whereas 


analysis, this change could be explained if the amount of negative feedback depended on the 
statistical properties of the input signal. However, given the different long-range behavior of 
spatial correlations for different ensembles, this would require to extend the model in order 
to include some spatial dependence of the feedback signal, for example, associated to the 
spiking activity of neurons that respond to stimuli presented in a shifted position. 

In the purely linear case, positive feedback always produces unstable dynamics. If an 
ON cell is subject to positive feedback, any positive stimulus fluctuation, no matter how 
small, feeds a reverberating loop where activity eventually diverges. Any negative stimulus 
fluctuation, in turn, eventually extinguishes firing altogether. An OFF cell is also unstable, 
with the opposite effect of positive and negative stimulus fluctuations. Such unstable systems 
are all-or-none (divergence or extinction) and have therefore not been studied here. 

An important finding of this paper is that negative feedback, which is usually assumed to 
exert a regularizing effect, can also produce unstable dynamics in the linear case. Negative 
feedback always diminishes the responses to slow stimuli, enhancing the band-pass charac¬ 
teristics of the filtering process. However, in certain conditions, these characteristics can be 
magnified dramatically, up to the point that the effective filter hr{uj) be sharply peaked at 
a specific frequency ui 0 . In the time domain, hr > {r) exhibits pronounced oscillations at this 
particular frequency. Feedback, hence, tr ansformed an intrinsic rec eptive held that acted as 
an integrator into an effective resonator Mat o and Samen go 200 8]. 

Neurons often display resonant properties. A widely accepted view states that oscillatory 
prope rties can eithe r stem from intrinsic cellular cha racteristics or from network interac¬ 


tions 


Buzsaki 2006 


Prescott et a 


two types of mechanisms coexist 


. 2008 


Wang 20iq. The important ingredient is that 


Hutcheon and Yarom 2000]: those attenuating high fre¬ 


quencies (typically, leak currents), and those attenuating low frequencies. Many processes 
can be invoked to attenuate low frequencies, both at the single-cell (e. g., leak filtering), 
and the network level (e. g., synaptic filtering). Previous studies of resonant behavior 
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have attempted to underst and resonances i n terms of such properties, invoking particu- 


lar subthreshold processes Izhikevich 2007, 

Gutfreund et al. 1995, 

Hutcheon et al. 1996a. 

Hutcheon et al. 1996b. 

Richardson et al. 200 

31, or specihc network interactions Wang 2010, 

Lampl and Yarorn 199" 

L Kondgen et al. 200 

3, Buzsaki and Wang 2012 . In addition, to 

avoid falling into quies 

ment with previous s 

icence, amplifying processes are sometimes also invoked. In agree- 

tudies on properties of neurons with spike-frequency adaptation 

Benda and Herz 2003, 

Gigante et al. 2007, 1 

3enda and Hennig 2008 

. Benda et al. 20101]. in 


this paper we have shown that negative feedback suffices to attenuate low frequencies. We 
make no assumptions about subthreshold properties, inasmuch as they produce an intrinsic 
receptive held that contains a maximal cutoff frequency. In addition, no ad-hoc mechanisms 
are required to produce amplification, as long as the baseline bring rate remains positive, 
for which a strong stimulus baseline so, or a strong spontaneous rate ho, suffice. Moreover, 
the formalism proposed here is general enough as to be equally applicable to network- 
mediated feedback currents, or to intrinsic adaptation currents. Our idealized approach, 
hence, provides a unified description of the mechanisms through which spike-triggered neg¬ 
ative feedback induces resonances. 

In order to test the validity of the idealized approach provided by perfect feedback, 
we also ran numeric simulations where a finite number N of individual feedback signals 
Xi{t) were triggered by actual spikes. Importantly, we concluded that the effective filtering 
characteristics do not depend on N. Perfect feedback, hence, can be safely used to study 
both self-inhibition (N = 1) and network-induced regulation (N > 1). The analytic study of 
the fluctuations in the feedback signal and in the response, however, shows that the amount 
of noise in the output spectrum diminishes as 1/N (Eqs. (122 ) ancl (l25l) h Moreover, the 
shape of the noise spectrum is given by a Cauchy distribution Johnson et al. 1994], 

The above conclusions hold for linear Poisson neuron models. Only for the linear case 
can the analytic approach be developed. The concept of receptive held, however, can also 
be extended to the case of linear-nonlinear Poisson models, or generalized nonlinear models. 
Although, strictly speaking, the presence of a nonlinearity does not allow us to employ linear 
methods, in the last section of this paper we demonstrated that by linearizing one of the two 
equations governing the system, an approximate description of the nonlinear case is possible. 
We showed that feedback still produces an effective receptive held that is narrower than the 
intrinsic one, and that resonances may also appear. Moreover, the nonlinear description is 
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also useful to understand how the divergences obtained in the purely linear case saturate at 
a finite value in the nonlinear description. 
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Appendix A: Appendix 

The convention used here to operate with the Fourier transform is 

r+oo 

I fit) e~ i0Jt df. (Al) 

— OO 

From this definition, the following properties follow: 

- The Fourier transform of a constant signal of magnitude r 0 is f(u) = \/27r r 0 <5(u;). 

- If a signal /(f) is equal to the convolution of two other signals g(t ) and h(t), then 
/M = g(u) h(u). 

- If a signal /(f) is equal to another signal but delayed, /(f) = g(t — A), then f(u>) = 
e~ iuA g(u). 

- If a signal /(f) is the derivative of another signal /(f) = dg/dt, then f{u) = iug{u). 

The Laplace transform, in turn, is defined as 

poo 

f(s ) = / /(f) df, (A2) 

Jo 

where s = cr + ico is a complex number. In the particular case where s is evaluated at 
a purely imaginary number (cr = 0), the Laplace transform is proportional to the Fourier 
transform for temporally positive functions. The Laplace transform, hence, can be seen as 
a generalization of the Fourier transform to the whole complex plane. Related properties 
hold: 
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- The Laplace transform of a constant signal of magnitude ro is f(s) = r^/s. 

- If a signal f(t ) is equal to the convolution of two other signals g(t ) and h(t), then 
f(s) = g(s ) h(s). 

- If a signal f(t ) is the derivative of another signal f(t) = dg/dt, then f(s) = s g(s) — 
9(0). 

Using these properties, it is easy to see that if the signals r(t) and x(t) are governed by 
Eqs. (0J) and (|8|) , the Laplace transform of r(f), for an initial zero feedback contribution 
x(0) = 0, is 


r(s) = 


According to control theory 


r 0 /s + h(s) si(s) 

+ !/ r d) 


(A3) 


Franklin et ah 1994], this feedback system becomes unstable 


when at least one pole of r(s) has positive real part. It is important to search for instabilities 
in the Laplace representation, since they may not be evident in the Fourier space. In addition 
to the fixed pole s = 0 given by the constant signal term, the poles of r(s) are those complex 
points s where the denominator of Eq. (1A3D vanish, that is, 


1 + 3 - 


h(s) 


= 0. (A4) 

S + 1/Td 

This equation holds in the complex plane, so both the real and the imaginary part of 
the equality must vanish. Feedback gives rise to unstable behavior whenever at least one 
solution of Eq. (1A4j) has positive real part. The onset of instability, hence, appears when 
the pole (or pair of conjugate poles) with largest real part crosses the imaginary axis, from 
left to right. At the crossing, cr = 0, so the real and imaginary parts of Eq. (1A4D become 


Im[/r(o;)] — Td to Re[h(o;)] = 0, (A5) 

1 + y/2n g r d Re[h,(o;)] = 0, (A6) 

where now, h(oj) is the Fourier transform of the intrinsic filter, obtained when evaluating 
the Laplace transform at a purely imaginary point and dividing by y/Zir. The frequencies 
satisfying Eq. (1A5I) are independent of the feedback strength, g. Note that, in addition, 
Eq. (1A6|) can only be fulfilled by those frequencies u>i, from the set of solutions to Eq. (1A5I) . 
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that result in Re [( 2 ( 0 ;*)] < 0. Even among these, if g is small, no frequency c o satishes 
Eq. (1 A6[l . and the system is stable (the pole with largest real part is on the stable semi¬ 
plane). For a critical value of feedback strength, the condition imposed by Eq. (1A6h can be 
finally reached and the system becomes unstable. This means that for feedback strengths 
beyond the critical value, there is at least one pole on the unstable semi-plane and therefore 
the inverse Fourier transform diverges. 


Appendix B: Appendix 


In this section, we describe the t emporal p rocess ing of slow stimuli. We arrive at the 


same results as the ones derived by 


Urdapilleta and Samengo 2009], but here we base the 


analysis on the spectral properties of filters. In the absence of feedback, 


"+OO 


r(t) = h 0 + 


h{r) s(t — r) dr, 


r( ui) = 


h 0 6(u) + h{uj) s(cu) 


(Bl) 


In the present context, a slow stimulus is one that does not contain high frequency 
components. In other words, the input/output relation of the cell is only determined by 
the lowest frequency components of the filter h(oj), since the higher frequencies are not 
explored. The stimulus, hence, has to remain fairly constant throughout the time scale of 
the filter (given by the non-zero portions in the temporal filters shown in Fig. [[]). If only 
low frequencies matter, we may take the limit oj —> 0. In this context, we prove that 


ON cells behave as low-pass filters, see Fig. [1]A.. For uj —> 0, the Fourier transform of 
these filters is 


h( u>) = 


-*+00 


yjJ — oo 
^ r +°o 


h(r) e~ iu}T dr 

h(r) [cos(cur) — i sin(a;r)] dr 


\ J —oo 

i r + °° i 

—f= / K r ) [l-iwr]dr = —= {H + iu H x ) 

v J —oo v 


(B2) 
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where H = h{r) dr and H\ — — f+°° r h{r) dr. The small angle approximation 
can be used again, 


h{uj) 


H 


H 


Hi 

1+iu n 


\p2/J . 


cos 


Hi . . . 

uj — | + 1 Sin 


UJ 


Hi 

H 


(B3) 

where we have defined the (positive) delay <5 0 = —H\ /H. From this expression, it is 
easy to check that the gain of the Fourier transform for ON cells at low frequencies 
is constant and equal to H/y/2n. In addition, the phase decreases linearly with u, 
starting at 0 with slope —So (in Fig. [TJAl , a linear graph instead of the semi-logarithmic 
one would clearly show this linear relationship). 



To relate h{uj) to the neural response in the temporal domain, we simply make use of 
the Fourier transform of a delayed signal (see Appendix |AJ , 


r(t) ^h 0 + J - 1 !H e~ iuS ° s(o;)} « h 0 + H s(t - S 0 ). (B4) 

According to Eq. (IBdjh the amplitude of the response is exclusively determined by 
H ; simultaneously, neural processing introduces a fixed delay 5 0 between the response 
and the stimulus. 


• ON biphasic cells behave as band-pass filters, see Fig. [Tp. In the limit u —> 0, h(u>) 
is still given by Eq. (ED- However, symmetric biphasic hlters satisfy H = 0, so to 
obtain a meaningful description we have to perform the expansion around uj — 0 up 
to the second order, 


" + 0O 


h{uf) 


V2n J 

1 c 


h{r) 


1 — i UJ T — -co 2 t 2 


1 

~2 


dr 


V2n \ 

i Hi uj 


H + % uj Hi - - uj 2 H 2 i , 


1 + i uj 


H 2 

2 Hi 


(B5) 
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where H 2 = r 2 h(r) dr. Proceeding as before and defining a corresponding 

(positive) delay <5i = —H 2 /(2Hi), we obtain 


h(oj) « 


e -/2 



+ i sin 




^ g * 71- / 2 g— 


k( i(7r/2 —oj5i) 


(B6) 


Hence, in the limit a; —>■ 0, the gain of symmetric biphasic hlters depends linearly on 
the angular frequency, with a positive slope Hi/\/2tt. In addition, the filter phase is 
a linearly decreasing function of c o, with intercept at the origin 7 t/ 2 and slope — Si. 
Both characteristics can be observed in Fig. [Ip. 

To relate h(u) to the neural response in the temporal domain, it is useful to recover 
the imaginary unit as a factor and rewrite Eq. (IB6I) as 


h(w) K -§= e-” 5 ‘ (i u). (B7) 

V Ztt 

In this case, given that the Fourier transform of a convolution becomes a product in 
Fourier space, the term (* co) can be effectively associated to the stimulus, s(o;). In 
turn, by applying the Fourier transform of a derivative, and proceeding as before, it 
is easy to check that r(t) is proportional to the delayed stimulus derivative, with a 
factor of proportionality given by Hi and a fixed delay 8\. Explicitly, 


r(t) « h 0 + T Hi e lujSl [i u: s(o;)] 


ho + Hi 


ds 

df 


(B8) 


J (t-Si) 


• Corresponding OFF cells behave analogously to monophasic and biphasic ON cells 
(see Figs. [Q3 and CEP). Specifically, the Fourier transforms of these hlters are exactly 
given by Eqs. (IB 2 li and (1B5I) . in the limit oj —> 0. The only difference with the previous 
ON cells is that the proportionality factors H and Hi for monophasic and biphasic 
cells, respectively, are now negative. The presence of a factor (—1) affecting the whole 
expression can be incorporated into a multiplicative term, exp(±27r), which shifts 
the phase of the Fourier transform in ±7r. This shift can be observed in Fig. Q] by 
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comparing phases of corresponding filters (monophasic ON/OFF filters and biphasic 
ON/OFF filters). Additionally, the factor (—1) does not affect the magnitude of the 
Fourier transform of ON or OFF cells. In the temporal domain, this factor simply 
means that the relationships obtained for ON cells are also valid, but associated to the 
n egati ve of the stimulus or its derivative. These results agree our previous analysis 
Urdapilleta and Samengo 2009 ]. 
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